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Abstract. We have developed a numerical code to study the deformation (e = (I zz — I xx )/I ZZl where la 
are the moments of inertia) of neutron stars in rapidly rotation in a fully general relativistic calculation. 
We have found that the deformation is larger, depending on the angular velocity, than is generally assumed 
for gravitational wave estimations. Calculations were performed by employing the Bethe-Johnson I EOS 
(equation of state) and a new set of models by the Frankfurt group including A hyperons for several choices 
of their coupling constants to ordinary nucleons. Possible implications for gravitational wave searches are 
briefly discussed. 

Key words: dense matter - equation of state - relativity - stars: neutron - stars: pulsars: general - stars: 
rotation 
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1. Introduction 

The problem of detection of gravitational waves is nowadays of great interest. There are many experiments 
in development to detect possible sources of gravitational waves. 

From the theoretical point of view several possible sources of gravitational waves have been advanced, like 
close binary systems, non-spherical collapsing stars, non-radially pulsating compact stars, wobbing pulsars, 
among others. We shall address in this paper, the particular case of the expected emission of gravitational 
waves from wobbling pulsars. 

The dimensionless amplitude of the gravitational waves, due to the wobbling neutron stars, depends 
on some poorly known parameters like the so-called "gravitational deformation" e = (I zz — I xx )/I zz and 
the wobble angle 9 W . Concerning the "gravitational deformation" e, values in the literature of the order of 
e ~ 10~ 3 — 10~ 4 are generally assumed by different authors, e.g. Pines & Shaham (1974), Zimmermann 
(1978), Shapiro & Teukolsky (1983), and Barone et al. (1988), to be adequate values for "slow" radio pulsars. 
It is worth mentioning that, even though implicit in the calculations, previous extensive studies on rapidly 
rotating neutron stars like those by Butterworth & Ipser (1976; BI76), Friedman et al. (1986; FIP86), 
Komatsu et al. (1989a, b; KEH89a, KEH89b), Eriguchi et al. (1994), Lattimcr et al. (1990; LPMY90), 
Cook et al. (1992, 1994a, 1994b; CST92, CST94a, CST94b), Bonazzola et al. (1993; BGSM93) or Salgado 
et al. (1994; SBGH94) have not addressed their results to the analysis of such a relevant parameter for the 
emission of gravitational waves. 

The aim of this paper is to show that the "gravitational deformation" e (see appendix for details) may 
be much larger than the values usually found in the literature, namely, e ~ 0.1 — 0.2 at rotation rates near 
the break-up value. These values are 2-3 orders of magnitude higher than those commonly adopted for the 
calculation of the gravitational emission rate of wobbling neutron stars. This is an important point since the 
actual pulsar population includes a large subpopulation for which e due to rotation would be found in the 
latter range, and this in turn enhances the prospects for detection of wobble radiation from the abundant 
galactic members, even assuming modest values for the wobble angle. 

It is important to note that, although we have performed a fully general relativistic calculation of e, the 
actual value of this parameter might be lower than we have found. This is due to the fact that an actual 
neutron star is much more complicated than the simple fluid that we have assumed and it should have, as 
it is generally accepted, a crust that behaves differently from a fluid under rotation. A more complicated 
question, however, has to do with the effective value of e once the star is precessing, because it is not clear 
how the structure of the star adjusts itself under this circumstances. The best way to calculate e in a fully 
general relativistic approach should take into account the crust and the precession of the star simultaneously. 
Such a model is obviously quite difficult to solve in a fully general relativistic approach. Therefore, we expect 
that the actual value of e may be lower than our present estimates, but our results indicate that it is probably 
higher than usually assumed. Since values of e <~ 10~ 3 — 10~ 4 may be more adequate for "slow" rotating 
sources, the presented results may be viewed as firm upper limits. 
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We remark that other works have already considered pulsars as possible sources of gravitational waves, 
among them Press & Thorne (1972), Zimmermann (1978, 1980), Zimmermann & Szedenits (1979), Thorne 
(1987), Barone et al. (1988), Nelson et al. (1990) and Finn & Shapiro (1990). In particular, some of the 
authors among others, Zimmermann (1978, 1980), Zimmermann & Szedenits (1979), Barone et al. (1988) 
and Nelson et al. (1990), addressed their studies to wobbling pulsars, not considering, however, a general 
relativistic calculation of e as we have done. We hope that our work can provide quantitative results for a 
better evaluation of the problem. 

Since the structure of a neutron star depends on the adopted equation os state (hereafter EOS) we 
should first specify the choices made on them. In the present work, we performed calculations using a 
new set of EOS's derived by Rufa et al. (1990), which include the contribution of hyperons in the strongly 
interacting matter. We have also calculated some models using the relatively stiff Bethe-Johnson I (hereafter 
B&JI) EOS in order to compare the accuracy of our numerical code with previous computations. 
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2. Numerical Code 

Our numerical code has been based on the approach developed by BI76, Butterworth (1976) and FIP86. 

The method developed by BI76 is a generalization of Stoekly's (1965) work on rotating Newtonian 
polytropes. 

To proceed, the distance element is first written in the form 

ds 2 = -e 2 »c 2 dt 2 + r 2 sin 2 dB 2 er 2 »(d^ ~ todt) 2 + e 2 ^ (dr 2 + r 2 d9 2 ) (1) 

with the metric functions u, B, u> and A independent of t and <j). 

Using the above metric in the Einstein equations we get differential equations involving v 7 B,uj and A. 
These equations were solved numerically using the Newton-Raphson method, for a given EOS and angular 
velocity ft. 

The Newton-Raphson technique requires, however, a known solution as a guess. We start the calculations 
from a spherical solution as a guess. This solution is easily obtained by solving the TOV equations for a 
given EOS and central density. As TOV equations are written in the Schwarzchild coordinates, we have to 
translate the results to the metric used here, which was done following Butterworth (1976). Thus, having 
the spherical solution as a guess, we take a small value of ft and obtain a model for this f2. Models for larger 
fi's are obtained using previous models, with lower values of SI. There is , however, an absolute maximum 
^ (Qmax) that a uniformly rotating star can support. This O mox is the Keplerian angular velocity (Clk), 
namely, the angular velocity of a particle in circular orbit at the equator. The maximum rotation rate, 
however, depends on the EOS and on the mass of the star (or central density). 

In fact, although the Keplerian angular velocity is an upper limit, Ipscr & Lindblom (1989a, b) showed 
that the maximum rotation angular velocity is actually below the Keplerian angular velocity. They argued 
that gravitational radiation instability implies that Slu m ~ (0.86 — 0.94)fi fe , where the viscosity, included 
in their calculations, works as a damping mechanism against instability, having therefore an important role 
on the value of Slu m . Other instability analyses were performed by Lindblom (1992), Lindblom & Mendell 
(1992) and Weber & Glcndcnning (1991), with similar results. 

For a given EOS and fixed (the injection energy of a unit mass particle lowered from infinity to the 
star), it is possible to get a sequence of models with < < fife. The choice of 0, for given EOS, defines 
the mass of the star. is related with the corresponding metric function in Eq.(l) by 

= e 2v ? (2) 

where v p is the value of v at the pole of the star. 

The structure equation is given by the equation of hydrostatic equilibrium which, for a uniformly rotating 
star, has a first integral 
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(3) 



where h(p) is the comoving entalpy density, v = (O — w)r sin#i3e _2ly is the velocity of the fluid relative to 
the zero angular momentum observer and c is the velocity of light. 

On the other hand h(p) can be obtained from the EOS, since we are dealing with zero-temperature 
matter, and the relation with pressure and energy density is given by 



Obtaining h(p), from the above equation, we can generate a table containing h,p and e. Once h(p) is 
determined from Eq.(3) we can get p and e by interpolation. We, in particular, have used a cubic spline 
interpolator. 

Thus its possible to get p(r, 9) and p(r, 9) determining, therefore, the star's structure for a given EOS, 
mass (M) and ft. 

FIP86 obtained several sequences for < O < ilk for different EOS's in the literature and fixed values of 
(3. However, by this procedure, for different values of O they obtained stars with different number of baryons. 
Since we were interested in following the variation of the star's structure as rotation increases, we have chosen 
to adjust ft in order to maintain fixed the total number of baryons along the sequence < O < fi^. 

In the method used here (Newton-Raphson) , the equations for v , B, and ui are linearized. For the 
resolution of these linearized equations, we transformed them into difference equations on a finite grid in the 
(r, fi) plane (where \i = cos 9) . The metric functions v, B and u) are even functions of /x and in this way the 
grid covers the interval < /z < 1. For the angular coordinate the spokes are taken from the Gauss-Legendre 
quadrature values \i\ = 0, ...it/. We particuraly choose I = 6. The Gauss-Legendre quadrature technique 
gives a good accuracy and , it is not necessary to include a large number of spokes. For the radial direction, 
we have used up to 120 spokes. In our calculations convergence to one part in 10 3 or 10 4 was required. 

Our numerical code was checked by calculating some models with B&J I EOS and the same star 
parameters used by FIP86. The comparison of the resulting properties is given in Table 1 and the agreement 
between both set of calculations is quite good, giving confidence to our numerical code. 



Hp) = Jc 



rp dp' 
e'+p' 



(4) 
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3. Results 



As we mentioned before, our models were based on the EOS's by Rufa et al. (1990). These EOS's represent 
a new approach concerning the inclusion of hyperons in the dense nuclear matter. Those calculations were 
performed in the relativistic mean field approximation and depend on two coupling constants of the theory, 
g^x and g a \ . These constants are not well determined experimentally and they were considered as parameters 
for the different set of EOS. In fact Rufa et al. (1990) found that is even possible to have self-bound A matter 
and, therefore, absolutely stable states at zero pressure. The hyperon species are expected to be present 
above the nuclear saturation density p and they may be very important to understand the cooling history 
of the stars since they may enhance the neutrino emissivity (see Pethick & Ravenhall 1992 and references 
therein). In Fig. 1 we show their EOS's calculated for different values of the coupling constants. For 
comparison, we plotted also the B&J I EOS. The relative importance of A's increases from EOS1 to EOS4 
and the EOS becomes softer for decreasing densities. Note that the EOS's are significantly different only 
for p > 6 x 10 14 gcm~ 3 . Table 2 gives the resultating non-rotating maximum neutron star masses for those 
EOS's as well as the respective values adopted for the vector and scalar coupling constants. We can see, as it 
would be expected, that the masses decrease as softer EOS's are used. Moreover, EOS4 allows a maximum 
mass of only 1.26M . This value is smaller than ~ 1.42M©, the masses of the neutron stars in the binary 
system PSR1913+16 and thus such low values for the scalar and vector coupling constants are excluded. 
Since EOS4 can be ruled out by these results, we have concentrated our efforts on the remaining ones. 

The properties of the stellar models are given in Table 3. Figure 2 shows the gravitational deformation 
£ as a function of angular rotation velocity for the calculated set of models, including also a B&J I EOS 
model. As mentioned earlier in the introduction, the values of e that we have obtained may be an upper 
limit, since an actual modelling should involve the detailed structure of the neutron star, in particular, its 
crust and the precession. 

Let us now estimate the value of the amplitude of the gravitational waves, h, for a wobbling pulsar 
considering the models that we have studied. The amplitude can be given by h = (I&ttGF/ c 3 ^ 2 ) 1 / 2 (see, 
e.g., Zimmermann 1978), where F — Lcw/^r 2 is the flux with Lew being the gravitational wave luminosity 
and r is the distance to the star. Law can be taken from, e.g., Zimmermann (1980) or Shapiro & Tcukolsky 
(1983) who provide us the gravitational wave luminosity for a wobbling star. Thus, we can write 



where: e is the gravitational deformation, 9 W is the wobble angle, I zz is the moment of inertia with respect 
to the rotation axis, Q, is the angular velocity in rad/s and r kpc is the distance to the star in kpc. It is 
worth mentioning that the above equation for h gives us only a rough idea (probably not better than the 
order-of-magnitude) of the true amplitude, due to approximations done in its derivation. The equation for h 
is, in fact, derived for a rigid newtonian object rotating free of external torques in the standard quadrupolc 




(5) 
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moment formalism, and it is not completely appropriate, therefore, for strong gravitational fields. 

In Table 4 we present the results of our calculations assuming 9 W — 10 ^, for the stars modeled in 
Table 3, at a distance rk pc = 20, rotating at the keplerian angular velocity. We have obtained h ~ 10~ 25 
which might be detected when, for example, the Caltech-MIT LIGO antenna (see, e.g., Thorne 1987 and 
Finn & Shapiro 1990) or the VIRGO antenna (see, e.g., Giazzoto 1989 and Bradaschia et al 1990) become 
operative but only if the wobble motion is continuously excited at least for several months. Otherwise, 
"spike" (impulsive bursts) with durations of <~ 1 s may result if the wobble is excited by some catastrophic 
phenomenon (e.g. phase transitions) but damps out on a quadrupole emission timescale (see, e.g., de Araujo 
ct al 1994). Similar values for h, but considering triaxial deformation of the pulsar, were obtained by Finn 
& Shapiro (1990). 
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4. Conclusions 

The present calculations of rotating neutron stars were performed using a new set of EOS's, including the 
presence of A's in the strongly interacting matter (Rufa et al. 1990). The resulting maximum masses for the 
neutron star configurations give bounds on the values of the scalar and vector coupling contants, suggesting 
9a\ ^ 0.43, otherwise too soft EOS's are obtained. 

If we take stars with the same mass for the Rufa et al. EOS's, it is not possible to distinguish them 
through their rotational behavior. In other words, the existence of fast millisecond pulsars is not at odds with 
exotic possibilities like self-bound A matter. The rotational properties of neutron stars were investigated for 
configurations having a constant number of baryons, condition required for obtaning consistent theoretical 
sequences which may be considered as an evolutionary path for decelerating objects. 

There is a potentially important application of this paper related to the fully general relativistic values 
derived for the gravitational deformations e, which may be two up three orders of magnitude (at best) 
larger than is usually entertained. We have shown that these values of e, at relatively small wobble angles, 
and observed rotation rates, might allow the detection of sources located in the Galaxy when the advanced 
generation of antennas becomes operative. 

Finally, it is worth mentioning that the calculations of e may be improved, since it is possible, in 
principle, to calculate it, in a full general relativistic approach, including in detail, e.g., the crust of the 
neutron star and the precession. In doing this we could have a lower value of e, as compared with the 
values that we have obtained in the present work. But these new values would still be higher than usually 
considered in the literature, because these last values are adequate for "slow" radio pulsars. 

It is our aim in the future to perform a more detailed calculation of e taking in account the crust of the 
neutron star and the precession with the use of our general relativity numerical code. 
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Appendix: The calculation of e 

The parameter e is defined as: 
e = (I zz - I xx )/I zz (Al) 
therefore, it is necessary to calculate the moments of inertia to obtain it. Before going into the details of 
our calculation of the moments of inertia it is important to characterize the gravitational wave amplitude h. 
Strictly speaking, due to the fact that the formula to calculate h is newtonian (or, at best, pos-newtonian) 
one might use, to be consistent, a newtonian (or pos-newtonian) calculation of the moments of inertia. 

Even if one had a definite formula for h in a fully general relativity version its calculation would be 
problematic, due to the fact that it is not possible to define I xx (or I yy ) invariantly, as done for the I zz 
through the definition: 

i zz = J/n (A2) 

(where J is the angular momentum and fi is the angular velocity). 

It should be kept in mind that the formula of h gives us only an order-of-magnitude estimate of the 
true amplitude, due to the various approximations in its derivation. It is worth mentioning that a newtonian 
calculation of the J^-'s should be, in principle, acceptable due to the fact that the equation used to calculate 
the gravitational wave amplitude is also derived for a rigid newtonian object. To calculate the ly's using 
the newtonian theory would be, on other hand, misleading, because we use fully general relativity theory to 
model our rapidly rotating stars. Thus, we decided to calculate all the Ty's through the formula: 
hi = J (r 2 S l3 - x lXj )dM (A3) 
where r is the radial coordinate of the metric used and the x^s are their projections, and dM is the 
gravitational mass element, namely 



= I Be 



2\-2v 



dM = (—T(j + Tf + T£ + T*)^gdr d6 d<j) (AA) 

2p+ (e+ P )(l + , 2 ) l +2rsindwB 2 e 2X- i A^ + P^\ r 2 sinedrded(f) 
1 — tr 1 — ir J 

where A, v, B and w are metric functions, and v, p and e are the velocity, the pressure and the energy density, 
respectively. 

This definition, is in fact, a newtonian- like way to calculate the moments of inertia, and therefore 
not invariantly defined, although it takes into account how the matter is distribuited throughout the star 
structure. 

We point out that our formula to calculate the Jj/s produces I zz ^ J/Cl, although, for the models that 
we have studied, the difference is < 10% throughout the whole sequence, i.e., < f2 < fife. 

Finally we note that, although the definition of the /^'s is not unique, the values for e should not 
change strongly if one adopts another definition, because, this quantity is an implicit monotonic function of 
the eccentricity after all. 
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Figure Captions 



Figure 1. Frankfurt and B&J I EOS's (see the text for details) 

Figure 2. Gravitational deformation £ as a function of SI for the models 1, 2, 3 
and B&J I. 
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Table 1. Comparison between our calculations and HP's calculations (in parenthesis). Where: f3 - injection 
energy, Q - angular velocity, p c - central density; M - gravitational mass, R eq - equatorial radius; u c /Sl - 
central dragging in units of angular velocity; T/W - ratio of rotational energy to gravitational energy, 
v eq /c - velocity of a comoving observer at the equator relative to a locally nonrotating observer, and 
e c = Wl — Rp/R?, q - eccentricity (where R p is the polar radius). 





CI 


Pc 


M/Mq 


Req 


LU c /Cl 


T/W 


Veq/C 






(raa/ s) 


(1U g 1 cm ) 




(km) 















1.00 


1.34 


12.3 












(0) 


(1.00) 


(1.32) 


(12.1) 












3000 


0.95 


1.33 


12.6 


0.43 


0.017 


0.14 


0.32 


0.671 


(3000) 


(0.95) 


(1.31) 


(12.5) 


(0.42) 


(0.017) 


(0.13) 


(0.27) 


(0.676) 


4030* 


0.91 


1.32 


13.3 


0.42 


0.033 


0.19 


0.49 




(4030) 


(0.91) 


(1.30) 


(13.3) 


(0.42) 


(0.034) 


(0.18) 


(0.47) 




5746^ 


0.79 


1.31 


16.2 


0.41 


0.087 


0.34 


0.73 




(5700) 


(0.77) 


(1.29) 


(16.9) 


(0.40) 


(0.093) 


(0.32) 


(0.74) 



* Angular velocity of the fatest pulsar known 
t Keplerian angular velocity 
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Table 2. Non-rotating maximum neutron star masses for the Frankfurt EOS's as well the respective values 
of the coupling constants g^x and g a \. 



EOS 


9<j\ 


9u\ 


M max /M 


1* 






2.17 


2 


0.80 


0.70 


1.98 


3 


0.43 


0.43 


1.54 


4 


0.0 


0.10 


1.26 



* No A's are present 
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Table 3. Results of our calculations for the Frankfurt EOS's and B&J I EOS. Where M Q is the baryon mass 
that is held constant ( within ±1%) through the calculations of each model. 



MODEL EOS M /Mq M/M q ft e e c I xx I zz R eq 

(rad/s) (10 45 5 .cm 2 ) (km) 



1, 2 1.06 1.01 3000 0.12 0.50 0.78 0.88 16.0 

/ and 3 1.08 1.02 3950+ 0.22 0.75 0.83 1.06 19.8 



1, 2 1.45 1.57 4030 0.15 0.52 1.01 1.19 15.8 

77 and 3 1.46 1.57 4900+ 0.24 0.75 1.10 1.45 19.2 



1 2.02 1.82 4030 0.10 0.43 1.01 1.12 14.4 

III and 2 2.02 1.84 5696+ 0.23 0.67 1.16 1.50 17.2 



1.59 1.45 4030 0.089 0.38 0.64 0.70 12.8 

BkJ I BkJ I 1.60 1.47 6200+ 0.25 0.75 0.75 0.99 16.6 



t Keplerian angular velocity 
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Table 4. Values of h using Eq. 5 for the models I-III and B&JI, for the Keplerian angular velocity, taking 
pulsars at T]~ pc = 20 and with W = 10~ 5 . 



MODEL 


h (xlCT 26 ) 


I 


7 


II 


17 


III 


22 


BkJ I 


19 



18 



